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Abstract 

We have studied the synchronization induced by periodic inputs applied to the 
finite iV-unit coupled bistable Langevin model which is subjected to cross-correlated 
additive and multiplicative noises. Effects on the synchronization of the system size 
(N), the coupling strength and the cross-correlation between additive and multi- 
plicative noise have been investigated with the use of the semi-analytical augmented 
moment method (AMM) which is the second-order moment approximation for lo- 
cal and global variables [H. Hasegawa, Phys. Rev. E 67 (2003) 041903]. A linear 
analysis of the stationary solution of AMM equations shows that the stability is im- 
proved (degraded) by positive (negative) couplings. Results of the nonlinear bistable 
Langevin model are compared to those of the linear Langevin model. 
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1 INTRODUCTION 



The Langevin model has been employed as a useful model for a wide range of stochastic 
phenomena in physics, biology and chemistry (for reviews, see [HE])- In the last decade, 
extensive studies have been made for the bistable Langevin model given by [3]-[20j 

dx- J N 

—t. = Xi - X f + -J2(x j -xJ+lW+ZiM+XiTliit), i = l,...,N, (1) 

where J denotes the diffusive coupling, I(t) an applied input, and £j(t) and Xifji(t) express 
additive and multiplicative noises, respectively [for detail, see Eqs. (jHJ)-((Z|)]- Refs. [3]-[T6] 
have studied the bistable Langevin model for a single element with N — 1 and J = in 
Eq. (JT]). Because it is difficult to study the coupled bistable Langevin model given by Eq. 
(TTJ with an arbitrary N, previous studies adopted the mean- field model given by [T7]-[2"U] 

— = (l-J)x-x 3 + Jx + I(t)+£(t)+xri(x), (2) 

which is valid for N = oo. Here x = (x) represents the time-dependent order parameter, 
and £(t) and xr](t) denote additive and multiplicative noises, respectively. The phase 
transition and synchronization in coupled bistable systems described by Eq. (j2J) have 
been investigated [T7]-[20j. 

A common approach to finite-N stochastic systems is to make direct simulations (DS) 
for the Langevin equation given by Eq. (TjQ) or for relevant Fokker-Planck equation (FPE) 
|21j . Despite a recent development in computers, it is not easy to perform DEs of Eq. 
(TTJ for appreciable values of N, which requires the computational time growing as iV 2 
with increasing N. For AT- unit Langevin equations, the FPE method leads to (N + 
l)-dimensional partial equations to be solved with proper boundary conditions, which 
is usually very difficult. As a useful semi- analytical method for stochastic equations, 
Rodriguez and Tuckwell [22] proposed the moment method in which the first and second 
moments of variables are taken into account. In this approach, original A^-dimensional 
Langevin equations are transformed to (A^/2)(A^ + 3)-dimensional deterministic equations. 
This figure becomes 65 and 5150 for AT = 10 and = 100, respectively. Based on 
a macroscopic point of view, Hasegawa [2"5l 12"^] has proposed the augmented moment 
method (AMM), in which the dynamics of coupled Langevin equations is described by 
a small number (three) of quantities for averages and fluctuations of local and global 
variables. The AMM has been successfully applied to studies on the dynamics of coupled 
stochastic systems described by the linear Langevin model [2"H I25J, FitzHugh-Nagumo 
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model [231 121], Ho dgkin- Huxley model [27] and networks [28]. The AMM is a semi- 
analytical theory which is useful for a study of the finite- iV stochastic system, whereas 
conventional mean-field- type approximations are applied only to the N = oo case |17]-|20j. 
Cortical processing, for example, is performed by many coupled populations of neurons 
and each population consists of finite numbers of neurons. Effects of cluster sizes on the 
information processing in brain may be clarified in the AMM [291 130] - 

Depending on the properties of elements forming a nonlinear system, we may clas- 
sify them into two types A and B. In the type A, elements are excitable units or self- 
oscillators such as FitzHugh-Nagumo and Hodgkin-Huxley models. In the type B, units 
that form the system cannot oscillate on their own unlike excitable systems. An example 
of the type B is a bistable element, which has been used in biology, chemistry and for 
a study of neural networks [21]- [S3]- Interesting and intrigue properties of the bistable 
Langevin model have been extensively investigated [3]-[2U]. Stationary probability distri- 
bution, first-passage time and the stochastic resonance for subthreshold periodic inputs 
in the bistable Langevin model subjected to cross-correlated additive and multiplicative 
white (or colored) noise have been studied. There are, however, still many unsolved 
basic problems. For example, dynamics of the finite-iV bistable Langevin model when 
time-dependent inputs are applied has not been well understood. The purpose of the 
present paper is to apply the AMM to the bistable Langevin model. We will investi- 
gate the dynamical response to external periodic inputs in the finite- iV bistable Langevin 
model subjected to additive and multiplicative noises with the use of the AMM. We apply 
suprathreshold inputs: deriving inputs are suprathreshold in the sense that they can in- 
duce the emergence of the transition of states from one of the bistable states to the other. 
We may investigate the synchronization of forced oscillations induced by suprathreshold 
periodical inputs in the type-B system, just as in the case of the type-A ensemble where 
the synchronization among self-oscillations are studied. We can examine effects on the 
synchronization of the size of systems, coupling strength and the cross-correlation between 
additive and multiplicative noises in a semi-analytical way within the AMM. The syn- 
chronization in the finite- N coupled bistable model is recently investigated by numerical 
methods |3"3] . 

The paper is organized as follows. In Sec. 2, the AMM is applied to coupled bistable 
Langevin model subjected to additive and multiplicative noise with the cross- correlation. 
Numerical model calculations are presented in Sec. 3. In Sec. 4, we make a linear stability 
analysis of the stationary solution. The dynamical properties of the bistable Langevin 
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model are compared to those of the linear Langevin model. Sec. 5 is devoted to our 
conclusion. 



with 



2 Augmented moment method 
2.1 Bistable Langevin model 

Generalizing the model given by Eq. (JT]), we adopt the iV-unit coupled bistable Langevin 
model given by 

^ = F{xt) + m + G( Xi ) th® + lt\t) + J«(t), (3) 

4 e) (*) = 4 EM*)"**®] (i = h-..,N). (4) 

Here F(x) = —dU(x)/dx, U(x) denotes the potential: G(x) is an arbitrary function of 
x: J expresses the diffusive coupling: Z (with Z = N — 1) stands for the coordination 
number: I^ e '(t) is an external input: and r]i(t) express zero-mean Gaussian white 

noises with correlations given by 

(Vi^Viit 1 )) = a 2 6 tJ 5(t-t'), (5) 

(6(*)6(0> = /? a M(*-0, (6) 
(ViWZiV)) = eapSySit-t), (7) 

where a and /3 denote the strengths of multiplicative and additive noises, respectively, 
and e the cross-correlation between additive and multiplicative noises. 

We will study the dynamical properties of the coupled Langevin model with the use 
of the AMM [231 121], in which the three quantities of /i, 7 and p are defined by 

fi(t) = (X(t)) = i$>(t)>, (8) 

i 

7(*) = ^E(M^)-Mi)] 2 ), (9) 
pit) = ([x(t)-Mt)] 2 )- (10) 

Here X(t) with = A^ -1 J2i x i(t) expresses a global variable, p its mean, and 7 and 
p denote fluctuations in local (xi) and global (X) variables, respectively. Equations of 
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motion for p, 7 and p are given by (the argument t being suppressed: for details, see 
appendix A) 

^ = fo + f2l + SU 1 2 + ^{a 2 [g g 1 + 3(g 1 g 2 + g 9 3 h} + ea^g 1 + 3g 3l )} + &\ 



dl „ , 9 , w 9 x 9 „ /2JiV 



(11) 

2/i7 + 6/37" + (0 + l)(2i 2 + 2g g 2 )a^ + 20ea^ 27 + (^^J (p - 7) + 

(12) 

^ = 2/ lP + 6/ 3 7P + (0+l)(^ + 2^2)« 2 p + 20ea^ 2 p+^, (13) 



with 



a 2 ^ + 2ea/% + 92l) + P, (14) 



where ft = (l/£\)d E F(p)/dx e , gi = (l/£\)d i G(p)/dx e , and = and 1 in the Ito and 
Stratonovich representations, respectively. The 0(7 2 )-order terms in Eqs. ffTTj) - ffT3]) are 
included with the use of the Gaussian approximation given by 



((^) 3 > ^ 0, (15) 

((^) 4 > ~ 3((5x t ) 2 ) 2 , (16) 

((6x t ) 2 (6 Xj ) 2 } ~ 3<(fe,) 2 ><(^) 2 >. (17) 

These terms play crucial roles for the bistable Langevin model although they are not 
necessary for the linear Langevin model [21]. Original N- dimensional stochastic equations 
given by Eqs. (EJ) and (jlj) are transformed to three-dimensional deterministic equations 
given by Eqs. (fTT]) - (jHj) . 

For the bistable Langevin model with 

F(x) = x-x 3 , (18) 

G{x) = x, (19) 

equations of motion are given by 

da o a 2 a eaB . . , . 

= ;U _ ;U 3_3 /i7 + _^ + _^ + /(t)) (20) 

^ = 2(l-3p 2 -3 7 )7 + 2a 2 7+(^)(p-7) + P, (21) 

^ = 2 (l-3p 2 -3 7 )p + 2a 2 p + ^, (22) 
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with 



P 



«V + 2ea(3p + (5 2 . 



(23) 



Eqs. (l20l) and (I2T!) with A^ = 1 (J = e = 0) are in agreement with results obtained for a 
single bistable Langevin model [H], while those with N = oo (p = 0) agree with results 
obtained for N = oo bistable Langevin model with a mean-field approximation |20j. For 
J = 0, Eqs. (J2TJ and (J22J lead to 

P - 1, (24) 

which expresses the central-limit theorem. 

A linear stability analysis of the stationary solution will be made in Sec. 4.1 with the 
use of the deterministic AMM equations given by Eqs. (I20|) - (l22|) . as in the case of the 
coupled FitzHugh-Nagumo model [26J. 

2.2 Synchrony 

In order to quantitatively study the emergence of a synchronized state of the ensembles 
defined by Eqs. ()3]) and (J4j, we first consider the quantity S'(t) given by 



When all variables are in the same state: Xi(t) = X(t) for all % (the completely synchronous 
state), we obtain S'(t) = in Eq. ( )25l) . On the contrary, in the asynchronous state where 
p = 7 /jV, it is given by S'(t) = 2(l-l/N)-f(t) = S' (t) [23J. We may define the normalized 
ratio for the synchrony given by [23] 



which is and 1 for completely asynchronous (S' = S' ) and synchronous states (S' = 0), 
respectively. 

We may alternatively interpret S(t) as the normalized mutual correlation given by 



(25) 




(26) 



S{t) 



at) 



(27) 



lit) 



with 



C(*) 



1 



J2T,lxi(t)-»(t)]Mt)-m] 



(28) 



NZ 




(29) 
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We note that S(t) = for J = [Eq. H}]. 
The equation of motion for S(t) is given by 




P\ 2 J 

-)S + ^r(ZS + l)(l-S), 



(30) 



with the use of Eqs. (|2Tj), (l22|) and ([26]) . 

3 Model Calculations 

We have performed model calculations, solving the AMM equations by the Runge-Kutta 
method with a time step of 0.01. Direct simulations for the iV-unit Langevin model have 
been performed by using the Box-Mueller algorithm and the Euler method with 
a time step of 0.0001. Results are averages of 1000 trials. 

Our model given by Eqs. (j3J) and (jlj) includes five parameters of N, J, a, (3 and e. 
We will investigate effects of N, J and e for fixed values of a = 0.1 and /3 = 0.1 otherwise 
noticed. 

Periodic pulse inputs 

We apply periodic pulse inputs at t > t\ given by 



with A = 1.0, ti = 50, T p = 100 and t w = 10 where Q(x) denotes the Heaviside function: 
Q(x) = 1 for x > and zero otherwise. Figure QJa), (b) and (c) show the time courses 
of /i(t), 7(i) and S(t), respectively, for iV = 10, J = 0.2 and e = 0.5. Initial values of 
Xi(t = 0) are set to be —1.0. When a positive pulse input I(t) shown by the chain curve 
in Fig. [T](a) is applied at t — 50 to the state which has been randomized by noise, the 
average value of \i is changed to about 1.0. When a negative pulse input is applied to the 
state with [i ~ 1.0 at t — 100, the state is switched back to fi ~ —1.0. In these switching 
process, the local fluctuation ^(t) and synchronization S(t) are transiently increased. In 
order to investigate the relation among /z, 7 and S, 7 and S calculated by the AMM are 
plotted as a function of \i in Fig. [T](d). It is shown that in the process of /i(t) changing 
from —1.0 to +1.0, S(t) has the maximum value at fi ~ 1.0. In the reversed process, 
S(t) has a maximum value at fi ~ —1.0. The maximum value of S(t) for a process from 



A Et 9 ^ - *i - kTpWi + kT P + t w -t) 



k 



eit-tt- fcr p /2)6(ti + kT p /2 + t w - t)}, 



(31) 
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H ~ —1.0 to /i ~ 1.0 is smaller than that for the inversed process, which is due to the 
introduced cross-correlation (e = 0.5) between additive and multiplicative noise. 

Similar fi — S plots for positive and negative correlations are depicted in Figs. [2](a) and 
(b), respectively. Figure [21^a) shows that with increasing e, the maximum value at fi ~ 1.0 
(yU ~ —1.0) is increased (decreased). The reversed behavior is realized for negative e, as 
shown in Fig. EJ^b): the fi — S plot for negative e is symmetric to that for positive e. 

Figures [3|a) and (b) show the J dependent fi — S plots for positive and negative J, 
respectively. From a comparison between Figs. |3](a) and (b), we note that the positive 
coupling is more effective than the negative one in increasing the synchrony: note that 
the vertical scale of Fig. |3]^b) is smaller than that of Fig. |3]^a). 

The maximum value of S, S max , is plotted as a function of N in Fig. @] which shows 
that the synchrony is more increased for larger J and smaller N. 

Sinusoidal inputs 

Next we apply sinusoidal inputs at t > ti given by 



with A = 1.0, t\ = 50 and T p = 100. Figures EJ^a), (b) and (c) show time courses of fi, 
7 and S, respectively, and Fig. [5^d) shows the relevant — S plot. From a comparison 
between Figs. [1] and [51 we note that the magnitudes of 7 and S for the sinusoidal input 
are about three times larger than that for pulse input. The \i — S plots for periodic pulse 
and sinusoidal inputs are similar besides their magnitudes. 

4 Discussion 

4.1 Stability analysis 

We will investigate the stability of the stationary solution of Eqs. ( 120]) - ( 122]) . from which 
the Jacobian matrix is given by 




(32) 



2 u 
1 — 3fi — 37 + — — 3/i 

-12/i 7 + 2(a 2 /i + ea/3) 2(1 - 3/i 2 - 67 + a 2 ) 
2 

12/i/9 + — (a 2 /i + ea(3) — 6p 







2JN 



2JN 



Z 



z 

2(1 - 3fi 2 - 37 + a 2 ) 



\ 



(33) 
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In the case of I = J = e = 0, we may analytically obtain stationary solutions and 
eigenvalues of the Jacobian matrix. Stationary solutions are divided into two cases A and 
B in which \i = and // ^ 0, respectively: the latter is further classified to cases Bi and 
B 2 as follows. 
(1) Case A 

Stationary solutions are given by 



7 








(34) 
(35) 



P 



(36) 



and relevant eigenvalues are given by 




(37) 



(38) 
(39) 



with 



D 1 = 



(l + a 2 ) 2 + 6/3 2 . 



(40) 



(2) Case B 1 

Stationary solutions are given by 




(42) 



(43) 



(41) 



and relevant eigenvalues are given by 




(44) 
(45) 



with 



D 2 = l-a 2 -y-3/3 2 . 
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(46) 



Table 1: Stationary solutions and eigenvalues within 0(a 2 ) and 0((3 2 ) 





oase A 


oase £?i 


oase r>2 







±(4-a 2 -3/3 2 ) 


|(a 2 + 3/3 2 ) 


7 


~(2 + 2a 2 + 3/3 2 ) 


i(a 2 + /3 2 ) 


^(4 + a 2 - 3/3 2 ) 


P 


^(2 + 2a 2 + 3/3 2 ) 


4^(« 2 + /^ 2 ) 


_J_ ( 4 + a 2_ 3/ 3 2) 


Ai 


-i(a 2 + 3/3 2 ) 


-2 + 2a 2 + 6fi 2 


a 2 + 3/3 2 


A 2 


-2(1 + a 2 + 3/3 2 ) 


-A -a 2 - 3/3 2 


-2 - 2a 2 - Q(3 2 


A 3 


-3£ 2 


-4 + 2a 2 + 3/3 2 


-3/3 2 



(3) Case B 2 

Stationary solutions are given by 

p 2 = ( 47 ) 

7 = ^(l + « 2 + yfj£), (48) 

P = ^( 1 + a2 + \/^)' («) 
and relevant eigenvalues are given by 

Ai, 2 = -2 + y/D~ 2 ± ^4-3/^2, (50) 
A 3 = -2 + a 2 + 2^. (51) 

When a and /3 are small, stationary solutions and eigenvalues are expressed in powers of 
a and 0: those in the lowest approximation are summarized in Table 1. We realize that 
the solution of A is stable while that of B 2 is unstable because its first eigenvalue of Ai 
becomes positive. The solution of B\ is stable if Re(Xi) < (see below). 

It has been shown that although results of the stability condition derived by the 
second-order moment method are not in good agreement with those of numerical method, 
it yields semi- quantitatively meaningful results [20]. Bearing this fact in mind, we will 
study the stability condition for the case of fi ~ ±1 (case Bi) with I — e — within the 
AMM by numerical methods. The solid curve in Fig. [6] shows the calculated boundary for 
J = 0.0 within which the stationary solution is stable (results for J ^ 0.0 in Fig. [6] will 
be explained shortly). Figure [7(a) shows the a dependence of the maximum eigenvalues, 
Xmax, for (3 = 0.0 and N = 10. For J = 0.0, A max becomes zero for a = 0.855 (/3 = 



10 



0.0), above which the stationary solution becomes unstable. When the coupling of J is 
introduced, the critical a value becomes 0.738, 0.968 and 1.106 for J = —0.2, 0.2 and 0.5, 
respectively. Figure [3(b) shows a similar plot of the (3 dependence of X max for a = 0.0 
and iV = 10. The solution becomes unstable for (3 > 0.577 with J = 0.0, a = 0.0 and 
iV = 10. With introducing J, the critical f3 value is changed to 0.518, 0.633 and 0.712 for 
J = —0.2, 0.2 and J = 0.5, respectively. Similar calculations of X ma x are made for a = 0.5 
with changing J. Squares, circles and triangles in Fig. [6] express results calculated for 
J = 0.5, 0.2 and -0.2, respectively, which are shown by curves for a guide of eye. Figure 
[6] shows that the stability of the stationary solution against additive and multiplicative 
noise is improved (degraded) by positive (negative) couplings. 

4.2 Effect of the symmetry of G(x) 

It is necessary to point out that the symmetry in the fi — S plot depends on the symmetry 
of G(x) for multiplicative noise. Indeed, in the case of G(x) = x which has the odd 
symmetry: G(x) = — G(—x), the asymmetry in the fi — S plot is obtained as shown 
in Figs. [2](a) and (b). However, if G(x) has the even symmetry: G(x) = G(—x), the 
asymmetry in the // — S plot is not realized. We have performed the AMM calculation 
for G(x) = x 2 — 1 with the even symmetry, for which equations of motion for //, 7 and p 
are given in appendix C. Calculated S(t) and j(t) are plotted against fj,(t) in Fig. [S^a) 
for e > and in Fig. Mb) for e < 0. With changing e, /i — S plot is little modified 
being symmetric independently of e, although magnitudes of the \i — 7 plot is changed. A 
comparison of Figs. M&) and (b) with Figs. [2(a) and (b) clearly shows that the symmetry 
of G[x) is important in studying the effect of the cross-correlation between additive and 
multiplicative noise. This fact may be applied to effects of the cross-correlation on the 
stationary probability distribution as shown below. From the FPE in Eq. ( lAlj) . the 
stationary distribution p(x) for I = J = in the Stratonovich representation (0 = 1) is 
expressed by 

N 

p(x) = Hpfa), (52) 

i=l 

with 

ln ^) - / WgW + wx) + W] dx - G) ln ^ G{x)2 + 2M) + ^ • 

(53) 
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Straightforward calculations for F( 3 , G(x) = x and G(x) 



x 



1 lead to 



p(x) oc (a 2 x 2 + 2ea(3x + (3 2 y a2+ ^ 1 - i ^ a4 - 1 / 2 

4ea(3x 2e(a 2 + /3 2 (3 - 4e 2 )) 



x exp 



a 2 x 2 



a' 



tan 



e(3 + ax 



p(x) 



OC 



[a 2 (x 2 - l) 2 + 2ta(5{x 2 - 1) + ^]-(i/2^+i/2) 

.tfe0 + a(x 2 - 1) % 



x cxp 



: tan 



for G(x) = x, (54) 



for G(x) = x 2 -l. (55) 



^Vl^ 2 V /^Vl^ 2 

Figures IHa) and (b) show the stationary distributions for G(x) = x and G(x) = x 2 — 1, 
respectively, with a = /3 = 0.5 for various e. For G(x) = x, an introduction of e yields the 
asymmetry in p(x), and p(a;) for a negative e is anti-symmetric with that for a positive 
e with respect to the x = axis. In contrast, for G(x) = x 2 — 1, p(a;) is symmetric 
independently of e, and an effect of a negative e is different from that of a positive e. The 
difference in G(x) reflects on various aspects of the bistable Langevin model such as the 
stationary distribution, the mean first-passage time and stochastic resonance, which have 
been conventionally calculated with the use of G(x) = x [3]-[T6]. 

4.3 The linear Langevin model 

It is worthwhile to compare the properties of the nonlinear bistable Lanvevin model to 
those of the linear Langevin model. For the linear Langevin model with F( 
(k: relaxation rate) and G(x) = x, we obtain equations of motion for fi, 7 and p in the 
Stratonovich representation as given by 



— = —KfJ, 



dfjL 

~dt 

dj 

~dt 
dp 

~dt 



a 2 11 eaB . , 



-2^7 + 2ct z 7 + 
-2np + 2a 2 p + 



(2JN\ 
, P 



){p-l) + P, 



N' 



with 



a 2 /r + 2ea(3p + (3 2 . 



Equations (I56 i) -( l59l) with e = agree with those previously obtained [2 
analytic expressions for the stationary state, given by 

21 - ea(3 



2k — a 2 



(56) 
(57) 
(58) 

(59) 

. We may obtain 
(60) 
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7 



P 



P 

2{k -a 2 - JN/Z) 
P 



1 + 



J 

Z(K-a 2 ) 



2N(K-a 2 ) 



yielding 



S 



J 



(61) 
(62) 

(63) 



J + Z{\-a 2 y 
where P is given by Eq. (159 j) with /i in Eq. (!6Up . 

Figure [TU1 (a) .(b) and (c) show time courses of /i, 7 and 5, respectively, for the linear 
Langevin model calculated with k = 1.0 and the same parameters for N, J, a, (3 and e 
as in Fig. [TJ The relevant fi — S plot is depicted in Fig. [T07 d). A comparison of Fig. 
[TUl with Fig. [1] shows that dynamical behavior of the bistable Langevin model is quite 
different from those of the linear counterpart. The nonlinearity in the bistable Langevin 
model plays an important role for the synchronization in the ensemble given by Eqs. ([3]) 
and OH). 

We may make a linear analysis of the stationary solution given by Eqs. (160]) -( 162]) . 
From Eqs. (I5"6"j) - (l5"8j) . we obtain the Jacobian matrix given by 



V 



a 



■k + 



2a 2 // 
2a 2 n 



N 



Eigenvalues are given by 



Ai 

A 2 
A3 





-2k + 2a 2 - 


k 4 



2JN 





2JN 
Z 

-2k + 2a 2 



\ 



(64) 



a 



2 ' 

-2k + 2a 2 
-2k + 2a z 



2JN 



(65) 

(66) 
(67) 



We note that eigenvalues in the linear Langevin model are independent of I and (3, which 
is different from those in the nonlinear bistable Langevin model. 



5 Conclusion 

We have studied the synchronization induced by periodic pulse and sinusoidal inputs in the 
iV-unit bistable Langevin model subjected to cross-correlated additive and multiplicative 
noise with the use of the semi-analytical AMM [231 El] • It has been shown that 
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(1) the synchrony is transiently increased when the mean value of state variables is 
switched from one stable state to the other which is induced by an external suprathreshold 
input, 

(2) the magnitude of synchrony is increased with increasing the coupling strength (J) 
and/or decreasing the system size (N), 

(3) The stability of the stationary solution against additive and/or multiplicative noise is 
improved by positive couplings but degraded by negative couplings, 

(4) the effect of the cross-correlation depends on its symmetry as well as a functional form 
of G{x) for the multiplicative noise, and 

(5) properties of the nonlinear bistable Langevin model are rather different from those of 
the linear Langevin model. 

The AMM [2"3"t |2~4"] may be applied not only to the type-A stochastic ensembles of excitable 
elements but also to the type-B ones consisting of nonexcitable elements. It is expected 
possible to apply the AMM to various types of coupled stochastic systems. In the AMM we 
may easily solve the three-dimensional deterministic equations of /x, 7 and p for periodic 
as well as non-periodic (transient) inputs, although its applicability is limited to the 
small-noise case which is inherent in the moment method. 
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Appendix A: Derivation of the AMM equations 

The Fokker-Planck equation for the Langevin model given by Eqs. (j3]) and (T4J) is given 

by m 

F(xi) + h + ^ [a 2 G f (x t )G(xi) + eafiG'(xi)] 



d 2 



p(x,t) J 

+ ^E^([« 2G (^) 2 + 2e «/ 3G (^) + / 32 ]^^)), (Al) 

where p(x,t) (x = {xk}), h = 1^ + I, G'(x) = dG(x)/dx, and = 1 and in the 
Stratonovich and Ito representations, respectively. 

With the use of Eqs. flAlj) . equations of motion are given by 



d -^L = (F(x l )) + (/ J ) + ^(G^G(x,)) + ^(G'(x 4 )), (A2) 
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d ^ tic j ) 



(xi F(xj)) + (xj F(xi)) + (xjj) + (xJi) 



y 

+ [a 2 (G(x,) 2 ) + 2ea/3(G'(x i ))+/3 2 ]% ) (A3) 

d(X ) 1 v v d(xi Xj) , . 

where X = X -1 Expanding Xj in Eqs. flA2p -f )A5|) around the average value of p as 

Xi = fi + Sxi, (A6) 

we obtain equations of motion for p, 7 and p given by Eq. ffTT]) - frH|) with the Gaussian 
approximation [Eqs. ( IT5]) - (TT7j) ]. 

Appendix B: The AMM equations for G(x) = x 2 — 1 
For the bistable Langevin model with 

F(x) = x-x 3 , (Bl) 
G{x) = x 2 -l, (B2) 



Eqs. (ITTl) - (fnj) yield equations of motion given by 

''" /i-/i 3 -3p7 + a 2 /i( / u 2 -l + 37) + ea/3 + /(t), (B3) 



dt 

^ = 2(l-3p 2 -3 7 )7 + 4a 2 7(3 / i 2 -l) + 2 e a/37+(^)(p-7) + J P, (B4) 

^ = 2(l-3p 2 -3 7 )p + 4a 2 p(3p 2 -l) + 2ea/3p+^, (B5) 
dt N 

with 

p = a 2 ( / , 2 _i)2 + 2ea/3(p 2 -l + 7)+/3 2 . (B6) 



15 



References 

[1] L. Gammaitoni, P. Hanggi, and F. Marchesoni, Rev. Mod. Phys. 70 (1998) 223. 

[2] B. Lindner, J. Garcfa-Ojalvo, A. Neiman, and L. Schimansky-Geher, Physics Report 
392 (2004) 321. 

[3] D. J. Wu, L. Cao, and S. Z. Ke, Phys. Rev. E 50, 2496 (1994). 

[4] C. Li, W. Da-jin and K. Sheng-zhi, Phys. Rev. E 52 (1995) 3228. 

[5] Y. Jia and J. R. Li, Phys. Rev. E 53, 5786 (1996). 

[6] Y. Jia and J. R. Li, Phys. Rev. Lett. 78, 994 (1997). 

[7] Y. Jia and J. R. Li, Physica A 252, 417 (1998). 

[8] L. Cao and D. J. Wu, Phys. Lett. A 260, 126 (1999). 

[9] D. C. Mei, G. Z. Xie, and L. Cao et al., Phys. E 59, 3880 (1999). 
[10] L. Cao and D. J. Wu, Phys. Rev. E 62, 7478 (2000). 
[11] J. Wang, L. Cao and D. J. Wu, Phys. Lett. A 308 (2003) 23. 
[12] C. W. Xie, D. C. Mei, and D. J. Wu, Eur. Phys. J B 33, 83 (2003). 
[13] Y. Jin, W. Xu and M. Xu, Chaos, Solitons and Fractals 26 (2005) 1183. 
[14] G-J. Zhang and J-X. Xu, Chaos, Solitons and Fractals 27 (2006) 1056. 
[15] Z. Ping, Chinese J. Phys. 44, 117 (2006). 
[16] H. Zhang, W. Xu and Y. Xu, Physica A 388 (2009) 781. 

[17] R. C. Desai and R. Zwanzig, Statistical mechanics of a nonlinear stochastic model J. 
Stat. Phys. 19 (1978) 1. 

[18] C. Van den Broeck, J. M. R. Parrondo, J. Armero, and A. Hern — 'andez-Machado, 
Phys. Rev. E 49 (1994) 2639. 

[19] M. Morillo, J. Gomez-Ordonez, and J. M. Casado. Phys. Rev. E 52 (1995) 316. 

[20] Y-M. Kang and Y-L. Jiang, Chaos, Solitons and Fractals 41 (2009) 1987. 

16 



[21] H. Risken: The Fokker-Planck Equation: Methods of Solution and Applications, 
Springer Series in Synergetics, Vol. 18 (Springer Verlag, Berlin, 1992) 



[22 
[23 
[24 
[25 
[26 
[27 
[28 
[29 
[30 
[31 
[32 
[33 
[34 
[35 
[36 

[37; 
[38 



R. Rodriguez and H. C. Tuckwell, Phys. Rev. E 54 (1996) 5585. 

H. Hasegawa, Phys. Rev E 67 (2003) 041903. 

H. Hasegawa, J. Phys. Soc. Jpn. 75 (2006) 033001. 

H. Hasegawa, Physica A 374 (2007) 585. 

H. Hasegawa, Physica D 237 (2008) 137. 

H. Hasegawa, Phys. Rev E 68 (2003) 041909. 

H. Hasegawa, Phys. Rev E 70 (2004) 066107. 

H. Hasegawa, Phys. Rev. E 75 (2007) 051904. 

H. Hasegawa, Physica A 388 (2009) 499. 

H. Sompolinsky, Phys. Rev. A 34 (1986) 257. 

M. Camperi, X. Wang, J. Comput. Neurosci. 5 (1998) 383. 

A. Koulakov, S. Raghavachari, A. Kepecs, J. Lisman, Nature Neurosci. 5 (2002) 775. 
J. M. Casado, J. Gomez- Ordonez and M. Morillo, larXiv: 0706 .15411 
S. Tanabe and K. Pakdaman, Phys. Rev. E 63 (2001) 31911. 

J. M. Sancho, M. S. Miguel, S. L. Katz, and J. D. Gunton, Phys. Rev. A 26 (1982) 
1589. 

R. F. Fox, I. R. Gatland, R. Roy and G. Vemuri, Phys. Rev. A 38 (1988) 5938. 
D. C. Mei, Y. L. Xiang and C. W. Xie, Phys. Scr. 74 (2006) 123. 



17 



Figure 1: (Color online) Time courses of (a) /x, (b) 7 and (c) S, and (d) 7 and S as a 
function of \i for the pulse input with N = 10, J = 0.2, a = 0.1, /3 = 0.1 and e = 0.5. 
Solid and dashed curves an (a)-(c) denote results of AMM and DS, respectively: the 
chain curve in (a) shows input I(t) shifted downward by - 2.0: arrows in (d) express the 
direction of time development. 



Figure 2: (Color online) e dependence of the fi — S plot: (a) for e = 0.0 (the solid curve), 
e = 0.5 (the dashed curve) and e = 0.9 (the chain curve), and (b) for e = 0.0 (the solid 
curve), e = —0.5 (the dashed curve) and e = —0.9 (the chain curve) with iV = 10 and 
J = 0.2. 



Figure 3: (Color online) J dependence of the fi — S plot: (a) The \i — S plot for J = 0.1 
(the dashed curve), J = 0.2 (the solid curve) and J = 0.5 (the chain curve), and (b) for 
J = —0.1 (the dashed curve), J = —0.2 (the solid curve) and J = —0.5 (the chain curve) 
with iV = 10 and e = 0.0. 



Figure 4: (Color online) The maximum value of synchrony S max against N for (J, e) = 
(0.5,0.0) (the solid curve), (0.2,0.0) (the dashed curve), and (0.2,0.5) (the chain curve) 
with a = 0.1 and j3 = 0.1. 



Figure 5: (Color online) Time courses of (a) /i(t), (b) j(t) and (c) S(t), and (d) 7 (the 
chain curve) and S (the solid curve) as a function of /i for the sinusoidal input with 
iV = 10, J = 0.2, a = 0.1, (3 = 0.1 and e = 0.5. Solid and dashed curves an (a)-(c) denote 
results of AMM and DS, respectively: the chain curve in (a) shows input I(t) shifted 
downward by - 2.0: arrows in (d) express the direction of time development. 



Figure 6: (Color online) The a-(3 phase boundary for J = 0.5 (the chain curve), 0.2 
(the dotted curve), 0.0 (the solid curve) and —0.2 (the dashed curve) calculated within 
the AMM with N = 10: the result for J = 0.0 is calculated by D 2 = [Eq. (j46j)], and 
marks (filled circles, squares and triangles) for J 7^ 0.0 show numerical results with curves 
plotted for a guide of eye (see text). 
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Figure 7: (Color online) (a) The maximum eigenvalue, X max , as a function of a with 
(3 = 0.0, and (b) X max as a function of (3 with a = 0.0 for J = —0.2 (dashed curves), 0.0 
(solid curves), 0.2 (dotted curves) and 0.5 (chain curves) with N = 10 (see text). 



Figure 8: (Color online) S and 7 versus /i for pulse input with G(x) = x 2 — 1, (a) S for 
e = 0.0 (the solid curve) and e = 0.5 (the dashed curve): 7 for e = 0.0 (the chain curve) 
and e = 0.5 (the dotted curve); (b) S for e = 0.0 (the solid curve) and e = —0.5 (the 
dashed curve): 7 for e = 0.0 (the chain curve) and e = —0.5 (the dotted curve) (N = 10 
and J = 0.2) (see text). 



Figure 9: (Color online) The stationary distribution p(x) for (a) G(x) = x and (b) 
G(x) = x 2 — 1 for various values of e with a — f3 — 0.5, the ordinate of (a) being different 
from that of (b). 



Figure 10: (Color online) Time courses of (a) /i, (b) 7 and (c) S, and (d) the \x — S plot 
for the pulse input applied to the linear Langevin model with iV = 10, J = 0.2, A = 1.0, 
a = 0.1, f3 = 0.1 and e = 0.5. Solid and dashed curves an (a)-(c) denote results of AMM 
and DS, respectively: the chain curve in (a) shows input I(t) shifted downward by - 2.0: 
arrows in (d) express the direction of time development. 
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